CellAreaInteger Function

private function CellAreaInteger(gridIn, i, j) result(cellarea)

Uses

Description compute area (m2) of a cell of a grid as a function of latitude defined by the position of cell in local coordinate system (row, column). Input grid of type grid_integer

Reference:

Sivakholundu, K. M., Prabaharan, N. (1998). A program to compute the area of an irregular polygon on a spheroidal surface, Computers & Geosciences, 24(8), 823-826.

Arguments

Type IntentOptional Attributes Name
type(grid_integer), intent(in) :: gridIn
integer, intent(in) :: i

row and column of cell

integer, intent(in) :: j

row and column of cell

Return Value real(kind=float)


Variables

Type Visibility Attributes Name Initial
real(kind=float), public :: cellWidht
logical, public :: check
real(kind=float), public :: ecc

eccentricity of spheroid
real(kind=float), public :: midLatitude
real(kind=float), public :: re

semi-major axis of spheroid
real(kind=float), public :: x
real(kind=float), public :: y

Source Code

FUNCTION CellAreaInteger &
!
(gridIn, i, j) &
!
RESULT(cellarea)

USE GeoLib , ONLY: &
!Imported parameters:
GEODETIC

USE Units, ONLY: &
!Imported parameters:
degToRad, squareKilometer

!arguments with intent(in):
TYPE(grid_integer), INTENT (IN)   :: gridIn  
INTEGER, INTENT (IN) :: i,j !!row and column of cell 

!Local variables
REAL (KIND = float) :: midLatitude !mid-latitude of cell
REAL (KIND = float) :: cellWidht !cell width in angular terms
REAL (KIND = float) :: re !!semi-major axis of spheroid
REAL (KIND = float) :: ecc !!eccentricity of spheroid
REAL (KIND = float) :: cellarea
REAL (KIND = float) :: x, y
LOGICAL             :: check
!------------------------------end of declarations-----------------------------


SELECT CASE (gridIn % grid_mapping  % system)
  CASE (GEODETIC)
    
    CALL GetXY (i, j, gridIn, x, y, check)
    IF (.NOT. check) THEN
       CALL Catch ('error', 'GridOperations',  &
       'cell out of grid calculating cell area' )
    ENDIF   
    
    midLatitude = y * degToRad
    cellWidht = gridIn % cellsize * degToRad
    re = gridIn % grid_mapping % ellipsoid % a
    ecc = gridIn % grid_mapping % ellipsoid % e 
    cellarea = re**2. * (1.-ecc**2.) * SIN(cellWidht)**2. * COS(midLatitude) / &
                    (1.-ecc**2.*SIN(midLatitude)**2.)**2.
  CASE DEFAULT
    cellarea = gridIn % cellsize ** 2.
  END SELECT

END FUNCTION CellAreaInteger